Stochastic current switching in bistable resonant tunneUng systems 
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Current- voltage characteristics of resonant-tunneling structures often exhibit intrinsic bistabilities. 
In the bistable region of the I-V curve one of the two current states is metastable. The system 
switches from the metastable state to the stable one at a random moment in time. The mean 
switching time r depends exponentially on the bias measured from the boundary of the bistable 
region Vth- We find full expressions for r (including prefactors) as functions of bias, sample geometry, 
and in-plane conductivity. Our results take universal form upon appropriate renormalization of the 
threshold voltage Vtti- We also show that in large samples the switching initiates inside, at the edge, 
or at a corner of the sample depending on the parameters of the system. 
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I. INTRODUCTION 

Recent advances in experimental techniques have made 
possible the study of fast stochastic processes such as dy- 
namic current switching in resonant tunneling structures. 
The electron transport in these devices has attracted a 
lot of attention since the pioneering work of Tsu and 
Esaki.-'^ The interest was further stimulated by the dis- 
covery of the phenomenon of intrinsic bistability^^^ in 
double-barrier resonant tunneling structures (DBRTS). 
Other resonant tunneling structures, such as superlat- 
tices, are also known to show bistable behavior. Re- 
cent experiments^"^" established that in the bistable re- 
gion one of the current states is metastable, and the 
switching to the stable state was studied. Both the mean 
switching time and its distribution function were mea- 
sured.^" 

The existence of intrinsic bistability is well understood 
theoretically. ^^"^"^ It was shown^* that in a certain range 
of bias, Vtti < V < Vth, for every value of V the current 
can take two different values, see Fig. 1. If one increases 
the bias starting from any value below Vth, the current 
follows the upper branch of the I-V curve shown in Fig. 1 
until V reaches Vth, where the current switches to the 
lower branch. On the other hand, if one decreases the 
bias from the values greater than Vth, the current follows 
the lower branch and then switches to the upper branch 
at Vth- 

The bistability can be understood by considering the 
potential profile of the DBRTS schematically shown in 
Fig. 2. If the level Eq in the quantum well is below the 
bottom of the conduction band of the left lead, tunneling 
into the well is not possible, and the current through the 
heterostructure is zero. In this case the charge in the well 
Q — 0- However, if a non-zero charge Q is added to the 
well, the level Eq rises due to the charging effects and may 
become higher than the bottom of the conduction band 
of the left lead. Then, another steady state of current is 
possible. In this state the current into the well from the 
left lead is compensated by the current out of the well 



through the right barrier. Thus, it is possible to have 
two different current states at the same bias. (See, e.g., 
points M and S on the I-V curve. Fig. 1.) 

The electric current in the device fluctuates, because 
the electrons tunnel in and out of the well at random 
moments in time. The resulting shot noise of current 
through the heterostructure gives rise to the metastabil- 
ity of some current states. The two solid lines in Fig. 1 
correspond to the most probable values of current at a 
given bias. These two branches are stable, i.e., any small 
fluctuation of current near a solid line will decay with 
time, and the current will return to its value at the solid 
line. The dashed line between Vth and Vth corresponds to 
the unstable state. Here any deviation from the dashed 
line which raises or lowers the current will switch the 
system to the upper or lower stable current state, re- 
spectively. 

Qualitative understanding of the metastability can be 
achieved by considering the system at a bias near Vth, 
e.g., point M on the upper branch of the I-V curve. 
Fig. 1. Then, as one can see from Fig. 1, a relatively 
small fluctuation can shift the current below the dashed 
line corresponding to the unstable state. If that happens, 
the system switches to the lower branch. The opposite 




FIG. 1: The I-V curve of the DBRTS^The bistable region is 
present in the range of bias between Vth and Vth- The bold 
dashed line corresponds to the unstable current state. 
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FIG. 2: Schematic potential profile of the double-barrier reso- 
nant tunneling structure. The structure consists of a quantum 
well separated from two leads by tunneling barriers. The elec- 
trons with three-dimensional wavevectors k and p fill all the 
states up to the Fermi energies Ef in the left and right leads, 
respectively. In the quantum well the motion of electrons 
in the z-direction is quantized, and the electrons with two- 
dimensional transverse wavevectors q occupy all the states 
up to the Fermi energy Ep. The inset shows the potential 
profile at zero bias. 



process is much less probable, since the distance from the 
lower branch to the dashed Hne is much larger than that 
from the upper branch. Therefore, the lower branch is 
stable, and the system remains in that state. 

The dependence of the mean switching time r on the 
bias was addressed theoretically in Ref. 15. It was shown 
that near the threshold voltage Vth the logarithm of r 
behaves as 



Inr cx 



(1) 



Here L is the size of the sample, rg oc cr^^^{Vth ~ V)^^l^ 
is a characteristic length scale, and a is the in-plane con- 
ductivity. In small samples, L <^ tq, the current switches 
from the metastable state to the stable one simultane- 
ously over the entire area of the device. On the other 
hand, in large samples, L 3> rg, the switching is initiated 
in a small critical region of radius tq. After the switching 
has occurred in that region, it extends rapidly to the rest 
of the sample. 

In this paper we show that if the sample is large, the 
switching can initiate not only inside, but also at the edge 
of the device. The latter process tends to be more effi- 
cient, since the exponential in the respective expression 
for T is smaller than in the case of switching far from the 
edges of the sample (Sec. V). On the other hand, the 
switching at the edge can be initiated anywhere along 
the boundary of the device, and thus the prefactor of the 
switching rate 1 /r due to these processes is proportional 
to the perimeter L. Similarly, the prefactor of the 



rate of switching inside the device is proportional to the 
area L^, which makes these processes more efficient in 
larger samples. 

We obtain analytically the full expressions for t, in- 
cluding the preexponential factors. Apart from the 
dependence on sample dimensions, the calculation of 
the prefactors reveals the non-trivial dependence of the 
threshold voltage Vth on the degree of disorder of the 
sample. Formal evaluation of the prefactors in the case of 
non-uniform electron density in the well results in ultra- 
violet divergences. Similar divergences appear in quan- 
tum field theory, where they are eliminated with the use 
of a renormalization procedure. The application of a 
similar technique to our problem leads to the renormal- 
ization of the threshold voltage which depends strongly 
on the conductivity of the quantum well (Sees. IV A 2, V) . 
Upon this renormalization In r in large samples acquires 
logarithmic corrections to its linear voltage dependence. 

The paper is organized as follows. In Sec. II we ob- 
tain the Fokker-Planck equation for tunneling in DBRTS 
which completely describes the electron transport in 
small samples. This equation enables us to find a simple 
result for the mean switching time in these samples. In 
Sec. Ill we derive the Fokker-Planck equation for the case 
of large samples which describes the dynamics of electron 
density in the well due to both the diffusion in the plane 
of the well and tunneling between the well and the leads. 
We use it to investigate the effect of weak density fluctu- 
ations on the decay of metastable current state in small 
samples (Sec. IV) and to study the switching in large 
samples (Sec. V). The application of our theory to the 
existing and future experiments is discussed in Sec. VI. 



II. FOKKER-PLANCK EQUATION FOR 
TUNNELING IN DBRTS 



The bistable current-voltage characteristic of DBRTS 
was studied theoretically in Refs. 11-14. The I-V curve 
shows the dependence of the average current on volt- 
age applied to the device. In addition, shot noise was 
studied in the regime of small fluctuations. ^^"^^ On the 
other hand, the switching between the branches of the 
I-V curve is caused by large fluctuations of current. In 
this section we use the model of Ref. 14 to derive the 
Fokker-Planck equation for tunneling in DBRTS, which 
accounts for these large fluctuations, and thus describes 
the switching. 

The model is illustrated in Fig. 2. The well is extended 
in x-y plane. The motion in z-direction in the well is 
quantized, and the well is assumed to have only one res- 
onant level of energy Eq. The two-dimensional wavevec- 
tors in the well are denoted by q. The left and right 
leads are three-dimensional; the wavevectors of electrons 
are denoted by k and p, respectively. The conduction 
bands in the leads are occupied up to the Fermi energy 
Ep. In typical devices Eq is of the order of Ep; for def- 
initeness we assume Eq > Ep. The temperature T is 
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assumed to be small compared to Ep and eV. The well 
is separated from the leads by two tunneling barriers with 
the transmission coefficients much smaller than unity. 

In Ref. 14 the tunneling through the double barrier was 
described quantum-mcchanically using the Breit-Wigner 
formula. The level widths with respect to the decay to 
the right and left leads Tl, there were eventually 
taken to be much smaller than all other relevant energy 
scales. We make this assumption from the beginning, 
and describe the electron transport through the barriers 
using the sequential tunneling approach. This method is 
an alternative to the use of the Brcit-Wigncr formula; it 
enables us to discuss both the I-V characteristic and the 
large fluctuations of current. 

In order to have a steady state of non-zero current 
in the device, the electrochemical potential in the well 
should lie between those in the left and right leads, i.e.. 



possible tunneling processes as 



= 2nJ2\tkfS{eV + E{k,)-Eo-e<t>), (6a) 



Tr = 2n^\tpf6{Eo+e<j)-E{p,)). 



(6b) 



To find Jl we use the Kronecker delta to remove 
the sum over k|| in Eq. (3), while the value of fcf = 
(2m/fi,^)(i?o + e(j) ~ eV) is fixed by the delta fimction. 
At k = (q, kz) and T the sum over q of /k(l — /q) 
can be easily evaluated, and gives {S/'in){kp — k^ — Qp) 
under the condition (2), where S is the area of the sam- 
ple. Then, the expression (3) can be simplified as follows 



Jl = 



h 



Sm 



{Ep + eV -Eo-e(l))-N 



(7) 



eV + Ep > Eo + e(j) + E^ > Ep. 



(2) 



Here Ef. = h^qp/2m is the Fermi energy in the well with 
m being the effective mass. Then, in the limit of low 
temperature the inequalities (2) dictate that the tunnel- 
ing is possible only in one direction, namely, from left to 
right, Fig. 2. The probability to tunnel through a barrier 
is given by the Fermi golden rule. The rates of electron 
tunneling into the well Jl and out of the well Jr take 
the form: 



qk 



Jr = 



xS{eV + E{kz)-Eo-e<f>), 

A-TT 



(3) 

-E{p.))- (4) 



Here E{k) = h k /2m; /k and /q are the Fermi occu- 
pation numbers in the left lead and the quantum well, 
respectively. In Eq. (4) we used the fact that the Fermi 
occupation numbers in the right lead /p — at energies 
above Eq + ecj), because Eq > Ep. Expressions (3), (4) 
include an additional factor of 2, which accounts for elec- 
tron spins. The matrix elements tp^ {tk^) describe the 
transitions between the resonant level in the well and 
the state with 2-component of the wavevector pz {kz) in 
the right (left) lead. The conservation of the transverse 
momentum is taken into account by Kronecker deltas. 

To simplify the expression for the tunneling rate (4) we 
use (5qp|| to remove the sum over p||. The remaining sum 
over q of Fermi function /q gives exactly the number of 
electrons in the well with a given spin N/2. Then Eq. (4) 
reduces to 



Jr = ^N. 

h 



(5) 



Here Tr is the level width with respect to tunneling into 
the right lead. We define the level widths for the two 



Here we used the expression N — Sqp/2TT for the total 
number of electrons in the well. Note that at eV > e(j) + 
Eq the level width (6a) vanishes, and thus Jl = 0. 

In the sequential tunneling approximation the average 
number of electrons in the well can be determined from 
the condition Jl = Jr, 



Sm Tl 
'^Tl+Ti 



-{Ep + eV-Eo-e(j)). (8) 



One cannot directly obtain from Eq. (8), since the po- 
tential (j) depends on the number of electrons in the well. 

Considering the barriers as two capacitors, one finds from 
electrostatics the following expression for the electric po- 
tential of the well (Fig. 2), 



V eN_ 
2" 2C' 



(9) 



Here we assumed for simplicity that the capacitances of 
the left and right barriers are equal to each other, and 
denoted the capacitance of each barrier as C. 

One can obtain the current-voltage characteristic of 
the DBRTS by repeating the following steps of Ref. 14. 
First, one notices that the level widths are energy depen- 
dent: 



Tl = gL^/Eo{Eo-eV + e4>), 
Vr = gR^/Eo{Eo + e(f>), 



(10a) 
(10b) 



where gL,R are dimensionless constants. Since F^ and Tr 
depend on (j), they are also functions of A^. Therefore to 
find N one has to solve the pair of equations (8) and (9). 
The latter leads to an equation on A^, which has three 
solutions in the bistable region. One of the solutions 
corresponds to the average number of electrons on the 
unstable branch, while the other two correspond to A'^ 
on the lower (A^ = 0) and upper stable branches. Upon 
substitution of N into Eq. (5) one finds the dependence of 
the average current on bias, i.e., the bistable I-V curve, 
which is schematically shown in Fig. 1. 
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To account for the noise, we go one step further and 
write the master equation for the time evolution of the 
distribution function P{N, t) of the number of electrons 
in the well N . In terms of the tunneling rates (3) and (4), 
the master equation for P{N,t) takes the form 

^-P{N,t) = P{N - l,t)JL{N - 1) 



dt 



+P{N+l,t)jR{N+l) 
-P{N,t)[JL{N) + jR{N)]. 



(11) 



The first two terms in the right-hand side of Eq. (11) 
account for the processes which increase the probability 
to have N electrons in the well, while the last term cor- 
responds to the opposite processes. 

In this section we consider the samples of large in-plane 
conductivity where the density in the well is uniform. 
Therefore, in the steady state of non-zero current the 
total number of particles in the well is proportional to the 
area of the sample. The linear dimensions of the sample 
are assumed to be large compared to the Bohr radius in 
the semiconductor. Thus the total number of electrons 
in the well is large, N ^ 1, and one can expand Eq. (11) 
in 1/N. Keeping the terms up to the second order, the 
master equation reduces^^ to 

|^(^'*) = -^[AWnN,t)] + \^^[B{N)P{N,t)]. 

(12) 

Here A{N) = Jl{N) - Jr{N) and B{N) = Jl{N) + 
Jb{N). Equation (12) is known as the Fokker-Planck 
equation, and is widely used for the description of various 
stochastic processes, see, e.g., Refs. 19,20. 

The stationary solution of Eq. (12) can be easily ob- 
tained: 



U{N) 



_ const 

= BiNf 

(13) 

Ref. 19. The extrema of function U{N) are determined 
by the condition A = Jl — Jr = 0, which we used above 
to find the average current through the device. There- 
fore each extremum of U{N) corresponds to one of the 
branches of the I-V curve. Outside the bistable region 
the current I{V) is uniquely defined, and U{N) has a 
single minimum. In the bistable region U {N) has two 
minima and a maximum, which correspond to the lo- 
cally stable current branches and the unstable branch, 
respectively [Fig. 1]. 

From the definitions of A and i?, it is clear that their 
ratio is independent of the area of the sample S. Using 
the expression (13) and the fact that N oc S, one can see 
that U is linearly proportional to S. Thus U{N) is an 
extensive quantity, and its dependence on N and S has 
the general form U{N, S) = Su{n), where n = N/S is the 
electron density. Since the area of the sample is large, we 
have U ^ 1. Therefore, the distribution function Pq is 
peaked sharply near the global minimum of U{N). 

The experiments^^^'^ studying the switching between 
the branches of the I-V curve are set up as follows. One 



starts at V < Vth, Fig. 1, where only one value of current 
is possible. In this case U (N) has only one minimum, as 
shown schematically by the dash-dotted line in Fig. 3a. 
If we increase the bias up to some value slightly above 
Vth, the function U{N) will acquire a new minimum to 
the left of the old one, see the dashed line in Fig. 3a. 
This corresponds to the appearance of the lower current 
branch of the I-V curve. The new minimum is a local 
one, and the main peak of the distribution function is 
still centered at the old minimum. Thus, the system 
remains on the upper branch of the I-V curve. Further 
increasing V, we transform U{N) to the shape shown 
schematically by the solid line in Fig. 3a. Here the right 
minimum of U{N) is a local one, and if we leave the 
system in this state for a sufficiently long time, it will 
eventually switch to the left minimum. To switch from 
the local minimum to the global one, the system has to 
overcome the barrier of height Ub, Fig. 3a. From the 
form of the distribution function (13) it is clear that this 
process takes a long time t cx cxp{Ub). To perform the 
measurement of the switching time from the upper to the 
lower current branch, one increases the bias to the chosen 
value over a time short compared to r, and then waits 
until the system switches to the lower branch. 

At the threshold the right minimum of U{N) disap- 
pears, and Ub = 0. Since U{N) is an extensive quan- 
tity, Ub grows rapidly when we move from the threshold 
into the bistable region, and the switching time becomes 
very long. Therefore, in order for the switching to oc- 
cur within a reasonable timeframe, the system should be 
close to the threshold. 

As the voltage approaches its threshold value, the max- 
imum at A^max and the local minimum of U (N) at N^nin 
(Fig. 3) move closer to each other, and at the threshold 
they coincide. At this point one can define a threshold 
electron density nth = A^max/'S' = N^in/S. In the vicin- 
ity of n = nth and V = Vui the function u{n) can be 
approximated by a cubic polynomial, 

7, 



u{n) w -a{n-nth) + ^{n-nthf+uth 



a = a{Vth-V). 

(14) 

Here the constant Uth is the value of m at n = nth and 
V = Vth- 

To derive Eq. (14) microscopically, one has to consider 
A{N) and B{N) on the upper branch of the I-V curve in 

the vicinity of the threshold Vth- An analytical calcula- 
tion of I{y) is possible^^ if the dimensionless parameter 



(15) 



is small, A <C 1. Here c = C/S is the capacitance per 
unit area. In Appendix A we extend this approach to 
find A{N), B{N), and the coefficients of expansion (14) 
at A < 1. 

The expansion (14) can be justified for any A in the 
spirit of the Landau theory of second-order phase tran- 
sitions.^^ The potential u is expected to be an analytic 
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(a) \ 



imated by a constant, 



(16) 




(b) 




V = Vth- dV 



FIG. 3: (a) Generic behavior of U{N) at different values of 
bias. Outside the bistable region U(N) has one minimum 
(top curve). Inside the bistable region the function U(N) has 
two minima and a maximum, which correspond to the locally 
stable current branches and the unstable branch, respectively 
(middle and bottom curves), (b) The sketch of U{N) for the 
model of Fig. 1. Solid line corresponds to a bias slightly below 
Vth, whereas dashed line depicts U{N) for the bias slightly 
above Vth- 



function of n and V. Thus, u can be expanded in Taylor 
series near the threshold, with n — nui playing the role of 
the order parameter. Since the local minimum and the 
maximum of u coincide at the threshold, both the first 
and second derivatives of u{n) vanish at V = Vth and 
n = nth- Therefore, the expansion starts with the third- 
order term. The sign of 7 is not important; we choose 
7 > 0, which corresponds to the behavior of U near the 
right minimum as shown in Fig. 3a. At V ^ Vth, the 
linear a{n — nth) and quadratic \l3{n — UthY terms are 
also present. Since a = l3 = QatV — Vth, we expect 
oc (Vth — y), /3 oc {Vth — V)- We keep only the linear 
term in the expansion, because the second-order term is 
quadratic in small parameter n — nth-, and therefore is 
small compared to the linear one. In order for u to have 
a local minimum at < Vth, the coefficient a should be 
positive. 

Near the threshold the function B{N) can be approx- 



In the case of constant B the Fokker-Planck equation (12) 
has been studied in detail. In particular, the exact ex- 
pression for the mean switching time can be obtained 
including the prefactor (Ref. 19, Sec. XIII. 2). In our no- 
tations it reads 



47r 



Sv'^"(^min)|C/"(7V^ax)| 



exp(C/f,). (17) 



For the potential (14) one can easily find the barrier 
height, Ub = 45'a3/2/3^i/2. xhc prefactor of Eq. (17) 
can be also straightforwardly evaluated, and one obtains 
the following expression for the mean switching time. 



2-K 



: exp 



4S[aiVth-V)f/^ 
3 



7 



1/2 



(18) 



where b = B/S is independent of the area of the sam- 
ple. This result obviously agrees with Eq. (1) for small 
samples {L <C ro). 

Expansion (14) is quite generic, and similar theoretical 
results were found in many different areas of physics. ^^^^^ 
In particular, Eq. (12) is also used to describe the mo- 
tion of a Brownian particle in external potential, where 
N plays the role of the coordinate of the particle. There- 
fore, the logarithm of the mean escape time of the Brow- 
nian particle from a local minimum of potential is also 
expected to obey the 3/2-power law. Recently this be- 
havior of the escape time was confirmed experimentally 
for the optically trapped Brownian particle.^® 

The lower branch of the I-V curve corresponds to the 
situation where the level in the well Eq + e(j) is below the 
bottom eV of the conduction band in the left lead. In 
this case = and B — —A — Jj^. Consequently, as 
N ^Owe have U{N) = 2N, see Eq. (13). Since N can- 
not be negative, U reaches its minimum at the boundary 
A'^ = of the range of allowed values of N, where the 
derivative U'{N) ^ 0, Fig. 3b. The non-analyticity of 
U{N) near the left minimum does not affect the calcu- 
lation of the time of switching from upper to the lower 
branch. Indeed, at a bias slightly below Vth, the function 
U{N) is analytic near its maximum and the local min- 
imum (solid line in Fig. 3b), and the description of the 
switching from the upper to the lower branch in terms 
of Eqs. (14) and (17) is correct. However, the situation 
is different for the switching from the lower to the up- 
per branch of the I-V curve. To study this process, we 
decrease the bias to the value slightly above Vth- The 
function U {N) for this case is depicted schematically by 
the dashed line in Fig. 3b. Here U{N) is non-analytic 
at its local minimum, and therefore we cannot use ex- 
pressions (14) and (17) for the switching time. The non- 
analytic behavior of U{N) is a consequence of crudeness 
of our model, in which the current on the lower branch 
is exactly zero. On the other hand, the experimentally 
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measured I-V curves show non-zero current on the lower 
branch. Thus, in a more detailed model which accounts 
for this non-zero current, the minimum corresponding to 
the lower branch of the I-V curve will be reached at non- 
zero A''. The discussion based on Eqs. (14) and (17) will 
then be valid. 



III. FOKKER-PLANCK EQUATION FOR 
TRANSPORT IN DBRTS OF LARGE AREA 

In section II we studied the decay of a metastable state 
in DBRTS under the assumption that the electron den- 
sity in the quantum well is uniform. Then the switching 
time T given by Eq. (18) grows exponentially with the 
area of the sample. Since electrons can tunnel at any 
point of the quantum well, the tunneling process creates 
a non-uniform electron density. On the other hand, the 
diffusion of particles in the well leads to spreading of the 
charge across the sample. In small samples the spreading 
is fast, and the density becomes uniform. In samples of 
large area the electron density may change significantly 
before the charge spreads over the entire well. In this 
case the switching between the two branches of the I-V 
curve is initiated in a small part of the sample, and the 
switching time is not exponential in the area S. 

In this section we generalize the Fokker-Planck equa- 
tion (12) to the case of non-uniform density n(r), where 
r = (a;, y) is a position in the well. In subsequent sec- 
tions this equation will be used to study the decay of a 
metastable state in DBRTS of large area. 



A. Equation for distribution function of electron 
density in an isolated quantum well 

We begin by considering the simplest case of a quan- 
tum well not coupled to the leads. At finite temperature 
the electron density in the well fluctuates and can be de- 
scribed by a distribution fmiction P{ri{r),t}. Here we 
derive the Fokker-Planck equation for the distribution 
function of electron density due to the in-plane diffusion 
of electrons in the well. In section III B we add the tun- 
neling through the barriers and obtain the Fokker-Planck 
equation for DBRTS of large area. 

We consider density fiuctuations at length scales much 
greater than the inelastic mean free path. These density 
fluctuations are slow in comparison with the energy re- 
laxation time in the well. Therefore the system is in a 
local equilibrium, and the distribution of electrons at any 
point in the well is given by a Fermi function. Note that 
the chemical potential in this Fermi function is deter- 
mined by the electron density, and therefore varies from 
point to point following the dependence n{r). 

Let us choose a time interval At much smaller than 
the relaxation time for P{n{r),t} and large in compar- 
ison with the collision time, so that the motion of elec- 
trons can be treated as diffusive. Then one can write the 



following equation for the evolution of the distribution 
function, 

P{n{r),t -\- At} - P{n{r),t} = J j drirfrs 

P{n{r) +(5ni2(r),i}VKAt(ri,r2;n(r) -|-5ni2(r)) 



-P{n{r),t}WAt{r2,ri;n{r)) 



(19) 



Here (5rii2(r) = (^(r — ri) — (5(r— r2) is the correction to the 
density n(r) due to the displacement of one electron from 
point r2 to ri; the probability density WAt(ri, r2; n(r)) 
describes diffusion of an electron from a point ri in the 
quantum well to point r2 during the time interval At. 
Since the diffusion rate may depend on the electron den- 
sity, Wa* is a functional of n(r). 

Expanding the first term in the right-hand side of 
Eq. (19) up to the second order in 5ni2(r), one obtains 
the following equation: 



AP{n,t} 



S 



X [WAt (ri , r2 ; n) - WAt (r2 , ri ; n)] 
1/5 5 " 



+ 



2 \5n(ri) (5n(r2) 
+WAt{r2,rr,n)]\p{n,t}. 



[WAt(ri,r2;n) 

(20) 



The probability densities Wai to diffuse from ri to r2 
and back are not independent. 



WAt(ri,r2;n)e" 



= WAt(r2,ri;n)e' 



-M2/T 



(21) 



Here fii and fi2 are the electrochemical potentials at 
points ri and r2, respectively. For the case of elastic scat- 
tering by impurities considered in Ref. 15 expression (21) 
directly follows from Eq. (8) of Ref. 15. Generalization of 
Eq. (21) to arbitrary scattering mechanism is discussed 
in Appendix B. 

In order to find /i we need to account for the interac- 
tions between electrons. We limit ourselves to the charg- 
ing energy approximation; the electron exchange and cor- 
relation effects are neglected. Then at low temperatures 
T <C Ep, the values of the electrochemical potential are 
found by adding the electrostatic potential e^n/c to the 
Fermi energy, 



Mi,2 = — n(ri 2). 
c 



(22) 



Here the effective capacitance per unit area c is defined 
by e^/c = e^/c -I- 1/u, and 1/ is the density of states in 
the well. 

In short time At an electron can only diffuse over a 
short distance, so that lui — ^2] <^ T. Therefore us- 
ing Eq. (21), one can expand the expression in the curly 
brackets in the right-hand side of Eq. (20) to the leading 
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order in (/ii — lJi2)/T, and with the help of Eq. (22) obtain 



dri(ir2< — 



cT \5n{vi) 6n{T2) 



x[n(ri) - n{v2)\ + 



(5n(ri) (5n(r2) 



xWAt(ri,r2;n)P{n,t}. 



(23) 



To proceed further we need an expression for the tran- 
sition probability density W^t- This quantity is affected 
by all the relevant processes of electron scattering, such 
as clastic scattering of electrons by impurities, elcctron- 
phonon and electron-electron scattering. Instead of ac- 
counting for all these processes explicitly, we express 1Fa< 
in terms of in-plane conductivity cr, which can in princi- 
ple be measured experimentally. Assuming that electron 
motion is diffusive, we conclude that the average square 
of the distance traveled by an electron during a short 
time interval is proportional to At, i.e., 



WAt(ri,r2;n)|ri -ral^dra = GAt. 



(24) 



Here the constant G is proportional to the conductivity, 

G = 4Tcr/e^, see Appendix C. 

At small A.t the transition probability density M^a* can 
be expanded as 



W^At(ri,r2;n) = 6{ri 



r2) + -^V^6{ri 



r2) + . 



(25) 

The physical meaning of the first term in this expansion 
is that electron remains at its initial position ri at At = 
0. Thus the second term is needed to account for the 
electron diffusion. The coefficient in the second term is 
found by applying the expansion (25) to Eq. (24). 

Equation (23) can be simplified significantly using ex- 
pansion (25), and eventually takes the form 



d_ 
dt 



P{n,t} 



dr 



Sn 



c 



■TV 



5_ 
Sn 



P{n.t}. 



This is the Fokker-Planck equation for the evolution of 
the distribution function of electron density. The first 
term in Eq. (26) describes the spreading of the charge 
in the well, whereas the second term accounts for the 
thermal noise. 

It is instructive to substitute into Eq. (26) the equi- 
librium distribution function Po{^}- The latter has the 
Gibbs form e~^/"^, namely. 



Po{n} = exp 



e^n^(r) 
2g~ 



dr 



Here the energy per unit area e = e n /2c is chosen 
in a way that reproduces the electrochemical potential 

/i = de/dn in the form (22). It is easy to check that 
Po{n} satisfies the Fokker-Planck equation (26). 



Combined Fokker-Planck equation for tunneling 
and diffusion 



In this section wc obtain the combined Fokker-Planck 
equation which incorporates both the tunneling through 
the barriers and diffusion inside the well. We begin by 
generalizing the tunneling Fokker-Planck equation (12) 
to the case of non-uniform electron density. This is ac- 
complished by dividing the plane of the well into small 
pieces, so that the density is uniform within each piece. 
In the absence of in-plane diffusion, the distribution 
function of electron density in the entire plane is given 
by the product of distribution functions of its pieces, 
P = YijPji^j}- Applying Eq. (12) to each piece we 
obtain the following Fokker-Planck equation for the dis- 
tribution function of the entire quantum well, 



^P = y^ 

dt ^ ^T^- 



P. 



The functions A[Nj) and B{Nj) are extensive quan- 
tities, and it is convenient to rewrite them as A{Nj) = 
ASa{n) and B{Nj) = ASb{n), where AS* is the area 
of each piece. Replacing the sum with the integral over 
the area of the sample and d/dNj with the functional 
derivative 5/dn{rj), we find the continuous form of this 
equation: 



d_ 
dt 



P{n,t} = dr 



Sn 



1 s 

-a(n(r)) + -— 6(n(r)) 



P{n,t}. 

(27) 

Let us now take into account the in-plane diffusion 
of electrons, which was discussed in Sec. HI A. Because 
the tunneling and diffusion are independent processes, we 
can add the right-hand sides of Eqs. (27) and (26) and 
obtain the combined Fokker-Planck equation for DBRTS 
of large area: 



d_ 
di 



P{n,t} = dr 



Sn 



V^n- 

c 



-a(n) + l-^Kn) 
2 on 



a n (5 



Sn 



P{n,t}. (28) 



This equation generalizes Eq. (26) to the case of a quan- 
tum well coupled to the leads. 

In the vicinity of the threshold Vih the function b{n) 
can be approximated by a constant b = b{nth)- In ad- 
dition, one can substitute 2a/b = 2A/B = —u'{n), c.f. 
Eq. (13). At bias near Vth the function u{n) is given 
by the approximate expression (14), and Eq. (28) can be 
rewritten as 



a + ^{n- nthf 



-2r]V^n + — 
Sn 



P{n,t}, 



(29) 



where we defined rj = a jch. In Eq. (29) we omitted the 
term proportional to the temperature, since it is negli- 
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gible at low T. (The exact criterion is discussed in Ap- 
pendix D.) Thus from now on we study only the effect of 
the shot noise due to the tunneling of electrons at high 
bias eV » T, whereas the thermal noise is neglected. 

The stationary solution of Eq. (29) is found by setting 
the left-hand side to zero, 



Po{n} 
F{n} 



-F{n} 



dr 



-a{n - nth) + - ^thf + vC^nf 



(30) 



The functional F{n} has two contributions: the first two 
terms account for the tmuicling, and the remaining term 
is due to the in-plane diffusion. 

The functional F{n} is similar to the free energy 
in Ginzburg-Landau theory of phase transitions, with 
n — nth playing the role of the order parameter. For the 
case of uniform electron density in the well, Vn = 0, the 
functional F{n} coincides with Eq. (13). If the density 
is non- uniform, the gradient term (Vn)^ appears in the 
expansion in addition to the terms from Eq. (13). This 
gradient term suppresses large variations of the electron 
density. 



C. Dimensionless Fokker-Planck equation 

For the following discussion it is convenient to 
parametrize the electron density n(r) in terms of a di- 
mensionless function z{p) that vanishes at the minimum 

of u{n), 



(a7)V4 ■ 



cbJarj 




Here the density at the minimum n„iin = nth + v^/t 
can be easily found from Eq. (14). The Fokker-Planck 
equation (29) in terms of z{p) takes the form: 



dP{z,t} 
dt 



xP{z,t}, 



Sz 



-VpZ + z ■ 



z' + 



where 



Srja 



UoSz 
(32) 



(33) 



IV. DECAY OF THE METASTABLE STATE IN 
EXTENDED SAMPLES 

In Sec. II wc obtained the expression for the mean 
switching time in DBRTS under the assumption of uni- 
form electron density in the well. This assumption is 
valid only if the dimensions of the sample are small com- 
pared to the length scale rg of the density fluctuations, 
Eq. (31b). If the sample is large, the fluctuations of elec- 
tron density must be taken into account. 

In Sec. Ill we obtained the Fokker-Plank equation (29) 
which describes the time evolution of the distribution 
function of electron density. Unlike Eq. (12) for the case 
of uniform density, this equation has an infinite number 
of variables, since the density is difFc;rent at every point. 

The most general form of the multidimensional Fokker- 
Planck equation is 



dPjx, t) 

di 



£P(x,t), 



92 



dxidxj 



A,(x).(35) 



Assuming that the system has a metastable state, one 
can consider its domain of attraction fi. The domain 
boundary dVt is a scparatrix of the drift field K. The 
mean time of the first passage out of the domain has 
been found in Refs. 29,30. For the process described by 
Eq. (35) the mean switching time is obtained as doubled 
mean first-passage time^° and takes the form, 



2^d'^a;Po(x) 
E/,orf5iEA,(x)Po(x)^ 



(36) 



Here Pq is the stationary solution of Eq. (35). The form 
function /(x) is a stationary solution of the adjoint equa- 
tion. 



9/(x) 
dxi 



0. 



£t/(x,i) = ^ |^ir,(x) +^ A,(x)^j 

(37) 

In addition, /(x) is defined to vanish at the boundary 
517 and reach /(x) ~ 1 well inside O. 

In subsequent sections we use the expression (36) to 
find the mean time of current switching in double-barrier 
structures. 



The stationary solution Pq of Eq. (32) is given by 



Po{z} 



F = Uo dp 



2 z^ z^ 



(34) 

One can see that the characteristic value of the func- 
tional F is given by Uq, whereas the characteristic size 
To plays the role of a typical length scale of stochastic 
fiuctuations of electron density n(r). 



A. Mean switching time in small samples 

In samples with linear dimensions small compared with 
ro the density fluctuations are weak. In this section we 
study their effect on the mean switching time. We will 
show that even these weak fluctuations can result in sig- 
nificant change of r. 
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1. Evaluation of the mean switching time 

In order to bring the Fokker-Planck equation (32) to 
the form (35) we present z{p) as an expansion 



(38) 



where 4>i{p) are the normahzed eigenfunctions of the 
Laplace operator, — V^(/<i(p) = ei(j)i{p). In particu- 
lar ^o(p) = ro/^S and eo = 0. Since there is no 
current flowing through the boundaries of the sample, 
the eigenfunctions must satisfy the boundary conditions 
h ■ ^</'i(p)|bomi(ja,i.y ~ ^' ^^^^^^ n is a unit vector normal 
to the boundary. The xo-coordinate corresponds to the 
average electron density in the well, whereas the other co- 
ordinates describe small fluctuations of the density. The 
eigenvalues a are numbered in order of increasing mag- 
nitude; ei ^ Tq/S* > 1. 

To obtain the x-representation of the Fokker-Planck 
equation we substitute the expansion (38) into Eq. (32) 
and find 



d 

dxi 



{ti + l- 2(t>oXQ)xi + 



1 d 
Uo dxi 



dOO r\ 
1 2 sr^ ^ o I 



' dxi' 



(39) 



where 



iijk = jdp(t>i{p)<j)J{p)<j)k{p)■ 
T}\c stationary solution Pq = e"^ in terms of Xi can 
be found by substituting expression (38) into Eq. (34). 
Then the functional F takes the form, 



F{x} = C/o 



-^(e, + l-2</.oa;o)a;| + 



2(l)Qxl 



1=0 
oo 



i,j.k—l 



(40) 



One can easily verify that exp(— i^{x}) solves the Fokker- 
Planck equation with C given by Eq. (39). 

The stationary probability density Pq is sharply peaked 
at the minimum of the functional F, i.e., at z{p) — 
(x = 0). Therefore, keeping terms up to the second order 
in Xi in expansion (40), we can evaluate the integral in 
the numerator of Eq. (36) in Gaussian approximation: 



minimum of F to its global minimum. The switching 
process is dominated by the paths which go through the 
vicinity of the lowest saddle point separating the domains 
of attraction of metastable and stable states. The bound- 
ary of the domain f2 lies exactly at the saddle point and 
is orthogonal to the direction of the steepest descent. 

The integral in the denominator of Eq. (36) is domi- 
nated by the saddle point of F. The latter is found from 
the condition SF/Sz = 0. This equation has an obvious 
solution Zs{p) = 1. In x-representation it corresponds 
to Xo = 1/^0 and ajj = for i > 1. Expanding expres- 
sion (40) near this point up to the second order in Xi we 
approximate F near the saddle point by 



F{^} ~ Uo 



1 



1 \2 1 °° 



(42) 

In small samples ei > 1, and therefore F has only one 
unstable direction xq, whereas all other directions are 
stable. One can see from Eq. (42) that in this approxi- 
mation the boundary 5f2 is the plane xq = l/4>o- 

Since the boundary dfl is orthogonal to the xq- 
direction, the sum over i in the denominator of Eq. (36) 
reduces to a single term with i = 0. Comparing Eqs. (35) 
and (39) one finds that Dij = {b^J^/Uo)5ij. Noting that 
Dij is diagonal, the sum over j also reduces to the only 
term with j = 0. 

To find df/dxo one needs to solve Eq. (37). Noting 
that £0 = and using Eq. (39), we can write the adjoint 
equation (37) near the saddle point as 



UqIxq- 



1 



+ 



d 
dxo 



dxo 



= 0. 



Solving this equation, we obtain 



dxo 



■2Uo ^-1 

TT 



(43) 



(44) 



Here the prefactor was found using the fact that / = 1 
inside the domain fl (i.e., at Xq — oo) and / = at the 
domain boimdary xq = l/0o- 

Using Eqs. (42) and (44) we can evaluate the integral in 
the denominator of Eq. (36) in Gaussian approximation. 
Then dividing the numerator (41) by this integral, we find 
the following expression for the mean switching time. 



tTo. 



(45) 



Here r is the switching time (18) obtained without the 
inclusion of density fluctuations. The latter give rise to 
the renormalization factor 



=n 



1 



1 



(46) 



/oo oo oo 
l[dxiP{^} = ll 
■°° i=0 i=0 



27r 



Uo{e^ + l)^ 



(41) 



In a multidimensional case in order to switch from the 
metastable state the system has to pass from the local 



To estimate the product Tq we assume a rectangular 
geometry of the sample with length L and width w. Then 
the eigenvalues are given by 



2„2 



IT r, 



m 
1? 



(47) 
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where n, m are non-negative integers. 

In small samples ^ 1 and the expression for In Tq 
can be expanded as 

InTo ^ - E' r- = -i E' T;^^^^ (48) 



n,m=0 



n,m=0 



where the prime in the sum means that the term with 
n = m = is excluded. 

The infinite sum in Eq. (48) is logarithmically di- 
vergent. However, since the diffusion picture is only 
valid at distances greater than the mean free path the 
wavevectors of the density fluctuations cannot^^ exceed 
Therefore, we need to cut the sum off at m < L/l 
and n < w/l. 

At L ^ w ^ the sum (48) can be approximated by 
a two-dimensional integral and yields. 



S \/S 



(49) 



Note, that although in small samples the area S is small 
compared to Tq, the effect of density fluctuations may 
become significant at I ^ \fS. 

In the case of strip geometry, w <^ L, we separate the 
sum into two parts, with n = and n > 0. The first part 
gives the sum of which can be explicitly evaluated 
and results in a small contribution L^/6rQ <C 1 to InTo. 
In the second part we approximate the sum over m by 
the integral with an infinite upper limit. Then neglecting 
terms (L/ro)^, we obtain the sum of 1/n. Cutting off 
this sum as discussed above, we find 



Lw w 



(50) 



For simplicity, from now on we will consider samples 
with w <^ L <^ y/S. 



2. Renormalization of threshold voltage 

Using Eqs. (45), (49) and (18) we find the following 
expression for the mean switching time in small samples. 



27r 



r 



• : 



4 S[a{Vth 
[3 7 



y)]3/2 



1/2 



m 

27rr/ I 



(51) 



The second term in the exponential of Eq. (51) represents 
the correction (49) due to the density fluctuations. 

Let us consider the regime when the magnitude of this 
term is larger than unity, but still small compared to 
the first term in the exponential of Eq. (51). Then this 
correction can be interpreted as a shift of the threshold 
voltage in formula (18). Indeed, substituting Vth — > Vth + 
^Vthi with the shift 



5Vth = 



477 ar] I 



(52) 



into Eq. (18) and expanding it up to the first order in 
5Vth we reproduce the result (51). In experiments the 
threshold voltage Vth is not known a priori. If one treats 
it as a fitting parameter, Eqs. (51) and (18) are equivalent 
up to the first order in 5Vth- 

The last term in the exponential of Eq. (51) formally 
diverges at / ^ 0. Similar divergences have been stud- 
ied in quantum field theory in the problem of the decay 
of the false vacuum. According to Eq. (14), the 
shift (52) of the threshold voltage is equivalent to adding 
a linear in the order parameter term —a5Vth{n — nth) to 
the integrand of the fimctional (30) . This corresponds to 
the standard in quantum field theory method of renor- 
malization of action, Refs. 16,17,31-35. Such renormal- 
ization procedure removes all the divergences. 

The origin of the renormalization of the threshold volt- 
age can be understood as follows. The "action" F de- 
scribes the so called (f)^ field theory in two dimensions, 
where (/) = (n — nth) is a scalar field. An alternative 
approach to the renormalization of this scalar field the- 
ory is to integrate out the fast modes (pp corresponding 
to large wavevectors, while keeping only slow modes (j)s 
with small wavevectors in the action F. One can find that 
the averaging of gives the sum of inverse eigenvalues 
of Laplace operator identical to (48), so that the term 
^(j)^ after the integration over the fiuctuations of the fast 
modes gives rise to j{(l)%)4>s = —a6Vth4>s- Physically 
this renormalization corresponds to the averaging of the 
switching rate over fluctuations of the electron density n 
in the well with characteristic scales between the mean 
free path and the sample size. 

Due to the renormalization of the threshold voltage the 
parameter a is modified as a ^ a + aSVth- Therefore, 
the quantities which depend on a, such that r and To, 
are also renormalized. More precise expression for r is 
given by Eq. (18) upon substitution of the renormalized 
a into it. On the other hand, the small corrections to 
the prefactor of r due to the renormalization are more 
challenging to observe experimentally, and for compari- 
son with experiment they can be ignored. 



V. MEAN SWITCHING TIME IN LARGE 
SAMPLES 

So far we studied samples of small area S' <C Tq. We 
found that the switching occurs when the electron den- 
sity at the saddle point is uniform, because the diffusion 
processes are fast and they smooth out all density varia- 
tions. In large samples, S ^ rQ, the diffusion is slower, 
and the system can reach the critical density in a small 
part of the well. After the switching occurs in that part, 
the switching process spreads rapidly throughout the en- 
tire well. In this section we study the switching time due 
to these nucleation processes. 

To find the expression for the mean switching time r 
in large samples we need to obtain the minimum and the 
saddle points of the functional F in Eq. (34). They can 
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be found using the condition 6F/5z = 0, i.e., 

-V'^z + z-z'^ = 0. (53) 

The boundary conditions for Eq. (53) should account for 
the fact that there is no current flowing through the 
boundaries of the sample. Since the current is propor- 
tional to the density gradient Vn, according to Eq. (31a) 
these boundary conditions take the form 



n ■ 



boundary 



0, 



(54) 



where n is a unit vector normal to the boundary. The 
trivial solution z{p) = gives the minimum of the func- 
tional F{z}, while the saddle points can be found as non- 
trivial solutions Zs{p) of Eq. (53). 



A. Nucleation processes in very large samples 

Let us consider the switching in an infinite sample, 
S ^ oo. Due to the symmetry of the problem, the 
solutions of Eq. (53) should be azimuthally symmetric. 
Placing the origin of coordinate system at the center of 
switching region and writing Eq. (53) in polar coordi- 
nates, we find 

<(p) + ^<(p) - z,{p) + zl{p) = 0. (55) 

This equation should be solved with the boundary con- 
dition Zs{p) = at /9 ^ 00, since otherwise F{zs} oc S, 
and the switching time r oc e^^^'^ will be infinite at 
S ^ oo. One can show that this condition is consistent 
with Eq. (54), that is 2;^(oo) = 0. Indeed, Eq. (53) can 
be interpreted as a Schrodinger equation for a particle in 



potential — Zs, i.e., — (V +Zs 



Therefore, Za{p) 



has the meaning of an eigenfunction of a bound state; its 



asymptotic behavior at large distances is Zs 



'Isfp, 



so that 2:^(00) = 0. The non-trivial solution of Eq. (55) 
with the boundary condition described earlier can be ob- 
tained numerically. The result is shown in the inset of 
Fig. 4. 

The main exponential dependence of mean switching 
time Ti in an infinite sample is given by e^^^"'^. Substi- 
tuting the numerical result for Zs(p) into Eq. (34), one 
finds^^ 



Ti exp 



8C 



7 



where the numerical constant 

"(Vz.)2 .2 



C 



dp 



2 



7.751. 



(56) 



(57) 



Equation (56) is the counterpart of the result (18) de- 
rived for small samples, 5* <C Tq. Because of the depen- 
dence of ro on V , see Eqs. (31b) and (14), both types of 
behavior can be observed in a single device by tuning the 
bias. At the crossover, Tq ~ 5, the two results coincide. 
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FIG. 4: The sketch of the density profile at the saddle point 
Zsi^p) corresponding to the solution of Eq. (53) with the 
boundary conditions (54). The precise radial dependence 
Zb{p) obtained by solving Eq. (55) numerically is shown in 
the inset. 



The problem of stochastic current switching is similar 
to the problem of finding the probability of spontaneous 
decay of a metastable vacuum near a Peierls transition 
point in (1-1-1) dimensional scalar field theory. The lat- 
ter problem was solved in Ref. 35, and the exponential 
factor in the result for the decay time is analogous to 
Eq. (56). On the other hand, the prefactor of the de- 
cay time is essentially different from r*, since we study 
the shot noise described by classical Fokker-Planck equa- 
tion, while the false vacuum decay problem is inherently 
quantum mechanical. 



1. Evaluation of the prefactor 

In a finite sample the switching can occur anywhere in 

the well, hence the prefactor of the switching rate 
must be proportional to the area S. Thus, while Ti has 
a large exponential, Eq. (56), its prefactor r* is propor- 
tional to l/S and can be small in large samples. There- 
fore, to fully understand the switching one needs to find 

T*. 

The time evolution of distribution function P{z,t} 
in large samples is given by the Fokker-Planck equa- 
tion (32). To evaluate the prefactor of the mean switch- 
ing time we again use the expression (36). The proce- 
dure is similar to the one for small samples described in 
Sec. IV A. However, the integration in Eq. (32) is now 
over a large sample, and therefore the density at the sad- 
dle point becomes non-uniform. Fig. 4. This significantly 
complicates the evaluation of the prefactor t* . 

We evaluate both integrals in Eq. (36) in gaussian ap- 
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proximation. As in Sec. IV A the integral in the numera- 
tor of Eq. (36) is dominated by the minimum of F{z} and 
is given by the expression (41). The denominator of (36) 
is dominated by the saddle point. Presenting z{p) near 
the saddle point as z{p) — Zs{p) + v{p), we obtain the 
expansion of F{z} in the form, 

F{zs + i} = F{zs} + ^ Jdpv{p)[-Vl - 2zs + l]i{p). 

(58) 

It is convenient to evaluate the integral in Eq. (58) by 
expanding 



"nm^rnm 



(p), 



(59) 



where (j>nm{p) are the normalized solutions of the eigen- 
value problem 

[- - 2^s (p) + 1] 4>nm (P) = ~Xnm4>nm (p) • (60) 

The boundary conditions for this equation are given by 
Eq. (54). 

Equation (60) can be interpreted as a Schrodinger 
equation for a particle in the attractive potential —2zs 
with energy Xnm — 1- Since the potential is azimuthally 
symmetric, we can separate the variables as (pnmip) = 
QnmiQp)'^ rnif) ■ The solutions for the azimuthal part 
'i^mi'fi) are given by e^™"^. Below it will be conve- 
nient to use their real combinations, cos nup and sinmcp, 
and introduce the following notations: ^'0(95) = I/V^tt, 
*m(^) = {cosrrvp)/y/TT for m = 1,2, ... , and *„(</?) = 
(sin m<^)/-v/7r for m = —1, —2, 

Substituting the expansion (59) into Eq. (58) and using 
the orthonormality condition for the eigenfunctions, we 
find 



f = F{z.} 



2 ^nm-^nm' 
n=0 m= — cso 



A 



(61) 



The discussion leading to Eq. (61) did not rely on the 
assumption of large sample size. In the case of small 
samples Eq. (61) reproduces the expansion (42), if one 
identifies A = e — 1 . This relation is easily understood by 
noticing that in small samples the density at the saddle 
point is Zs{p) = 1. Comparing the definition of e given 
in the paragraph after Eq. (38) with Eq. (60), where we 
substitute Zs = 1, we reproduce A = e — 1. 

The form of Eq. (61) suggests that in the case of large 
samples it is more convenient to evaluate the integral in 
the denominator of Eq. (36) using variables Xnm rather 
than Xnm- Since the eigenfunctions (j)nm and are nor- 
malized, the expansion coefficients Xnm are related to co- 
efficients Xmyi of expansion (38) via an orthogonal trans- 
formation. The Jacobian of this transformation equals 
unity, and therefore the integration over Y\ dXnm in the 
denominator of Eq. (36) can be replaced by the integra- 
tion over Y[ dXnm ■ 



In order to evaluate the integral in the denominator 
of Eq. (36) in the x-representation, we need to find the 
eigenvalues Xnm of Eq. (60). All Xnm are positive with 
the exception of one negative eigenvalue, Aoo < 0, and 
two zero eigenvalues, Ao,i = Ao,-i = 0. Numerical solu- 
tion of equation (60) yields Aqo ~ —1.648. This negative 
eigenvalue is associated with unstable deviation from z^ 
corresponding to the motion over the saddle point. In 
Eq. (36) the boundary dfl of the domain of attraction 
of the metastable state is orthogonal to xoo-direction, so 
that the integration in the denominator is performed only 
over the positive and zero modes. Since each positive Xnm 
corresponds to a Gaussian integral, the integration over 
them is straightforward. The integration over the zero 
modes is more challenging; to perform it we first need to 
understand their physical meaning. 

The existence of two zero eigenvalues is due to the 
translational invariance of the functional F{z} with re- 
spect to any shift of the center of the switching region in 
the plane of the quantum well. The two zero eigenvalues 
correspond to two orthogonal to each other directions in 
the plane along which such a shift can be performed. In- 
deed, a small shift Ap of the center of switching region 
results in the following small change in the saddle point 
density, 

dz dz 
Zs [p + Ap) -Zs{p) = ^Ap^ + -^Apy. (62) 
^Px ^Py 

One can check by differentiating Eq. (53) with respect 
to px,y that the derivatives dzs/dpx,y are solutions of 
Eq. (60) with A = 0. Furthermore, dzs/dpx ~ z'^{p) cos ip 
and dzs/dpy = z'^{p)s\Tnp, so the azimuthal quantum 
numbers corresponding to zero modes are m = ±1 in 
our notations. Thus we conclude that dzs/dpx,y = 
co0o,±i(p), where cq is a constant. "^^ 

Substituting these expressions for dzs/dpx,y into 
Eq. (62) and comparing it with the expansion (59), we 
find that the coefficients corresponding to zero modes are 
a;o,±i = Co^Px,y Thus the integral over the zero modes 
^0,1 and xo,-i amounts to the integration over the possi- 
ble positions of the center of the switching region in the 
sample, 



jdxo,i jdxo,-i=cl j d{Apx) j d{Apy) = Q^. (63) 

Here the constant Cg was found using azimuthal symme- 
try of Zs and the fact that the eigenfunctions (^nm{p) are 
normalized, 

cl = J (Ij) ' dp=lj [Vzs{p)r dp = C. (64) 

The relation between the last integral and the constant 
C defined by Eq. (57) is proven in Appendix E. 

To find the denominator of Eq. (36) in the x- 
representation we also need Dij and df/dxoQ. They 
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can be obtained from the x-representation of the Fokker- 
Planck equation for large samples. Substituting z{p) = 
Zs{p) + v(p) with V in the form (59) into Eq. (32) and 
using the orthonormality condition for the eigenfunctions 
4>nm{p), we obtain the Fokker-Planck equation P = CP 
with 



E 



dx. 



(65) 



Here we neglected the terms of higher orders in Xnm- 
One can easily check that the solution of the stationary 
Fokker-Planck equation LPq = is Pq = with F 
given by the Gaussian approximation (61). 

Comparing Eqs. (35) and (65) we conclude that Dij = 
{by/^/Uo)Sij. To find df/dxoo we need to solve Eq. (37) 
with jC given by (65), that is 



-Aooioo + 



1 d 



Uo dxoQ 



df 
dx 



0. 



00 



(66) 



Solving it with the conditions / = 1 inside the domain fl 
(i.e., at xq — oo) and / = at the domain boundary 
xqo = 0, we find that at the saddle point df/dxoo = 
-(2|Aoo|C/o/7r)i/2. 

Substituting Eq. (41) for the numerator of Eq. (36), 
and Eqs. (61) and (63) along with the expressions for 
Dij and df/dxoo into the denominator of Eq. (36) we 
reproduce the result (56) with the prefactor given by 



(67) 



Here the product T'" excludes the factors corresponding 
to the three non-positive eigenvalues Xnm- The coeffi- 
cients Xnm denote the parameters l + c; used in Sec. IV A. 
They coincide with the eigenvalues of the Schrodinger 
equation (60) in the absence of the attractive potential 
-2zs{p). 

To evaluate the infinite product T'" we need to find the 
continuous spectrum of Eq. (60). The radial part of (f)nm 
oscillates as a function of p with the wavcvector (j„. The 
phase of these oscillations at p ^ oo is shifted by 6m{Qn) 
due to the scattering in the attractive potential —2zs. 
The cigcnvahics of the continuous spectrum are expressed 
in terms of these scattering phase shifts as follows 



1 + 9. 



/7rn\' 



i/R) 



7rn 



(68) 



This result is derived for a round sample of dimensionless 
radius i? » 1; the derivation and the expression for the 

phase shifts Sm arc given in Appendix F. The expression 
for the eigenvalues Xnm is given by Eq. (68) with 6m = 0. 

It is convenient to calculate the logarithm of T'", 
thereby transforming the product over n and m to a sum. 
Taking the large sample limit, R oo, we replace the 



sum over n by an integral over q = nn/R. Then expand- 
ing the integrand in small parameter Sm/n, we find 



InT"'~- 



qdq 

l + g2 



(69) 



To investigate the convergence of the integral we need 
to evaluate the sum of the phase shifts at large q. This is 
accomplished with the help of the following "Priedel sum 
rule" 



2C 



(70) 



proven in Appendix F. The asymptotic behavior (70) of 

the phase shifts implies that the integral in Eq. (69) di- 
verges logarithmically at g ^ oo. This ultraviolet diver- 
gence signals that T'" is determined by a large wavevec- 
tor cutoff or, equivalently, by some short distance scale. 
An analogous divergence appeared in the prefactor of the 
mean switching time in small samples. Sec. IV A. There 
we have shown that this short distance cutoff is of the 
order of the mean free path I. Following the same recipe, 
we cut off the integral in Eq. (69) at g ^ ro//, and with 
logarithmic accuracy find 



InT" 



¥-(t)- 



(71) 



Substituting this result into Eq. (67), we obtain the pre- 
factor 



1 



bSa^ \ro J 



I 



2C/7I 



(72) 



This expression completely describes the parametric de- 
pendence of the prefactor of the mean switching time in 
large samples. On the other hand, because of the ul- 
traviolet divergence of T'", the numerical coeflScient in 
T* cannot be determined without detailed treatment of 
charge transport at short length scales. 



2. Renormalization of threshold voltage in large samples 

Expression (72) for the prefactor r* implies that in 
large samples the switching rate diverges at I ^ 0. 
In Sec. IV A 2 we encountered the same problem while 
considering small samples. There it was shown that the 
dependence of t on the mean free path can be absorbed 
into the definition of the threshold voltage Vth- Following 
the same renormalization technique, one can shift the 
threshold voltage Vth by the amount 



th 



-J-^ln^, 
47r ar] I 



(73) 



chosen in such a way that the resulting correction to the 
exponential in Eq. (56) cancels T'" in the prefactor r*, 
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see Eqs. (67), (71). The renormalized result for the mean 
switching time then takes the form 



1 



bSa% 



cxp 



7 



aR = a + aSVth. (74) 



This expression is equivalent to Eqs. (56), (72) up to 
the correction in the prefactor n due to the substitution 
a aR. 

The characteristic length scale ro is sensitive to the 

position of the threshold vohagc, so its value has to be 
renormalized. Since the size of the critical nucleus and 
5Vth are coupled to each other, Eq. (73), they should be 
evaluated self-consistently: 



tr 
5Vth 



-^^in!:^. 

4ir arj I 



(75a) 
(75b) 



To find SVih one; can solve the system of equations (75) 
iteratively starting with tr = vq. The result (73) then 
should be understood as the first iteration of Eq. (75b). 

Upon the substitution of the shift (75b) into Eq. (74), 
the logarithm of the switching time Tj acquires an ad- 
ditional logarithmic dependence on voltage due to the 
bias-dependent rcnormalization of Vth- This dependence 
is physically meaningful and can, in principle, be tested 
experimentally. However, these corrections to the volt- 
age dependence (56) of r, are small, and to the leading 
order In Ti is still linear in voltage. 



B. Nucleation near sample boundaries 

In Sec. VA we studied the nucleation processes in 
very large samples assuming that the switching initiates 

far from the boundaries (i.e., at distances significantly 
greater than ro). In this section we show that the switch- 
ing can be more effective when it is initiated near the 
boundaries of the sample and evaluate the mean switch- 
ing time for such processes. 



1. Nucleation at a smooth edge 

To study the nucleation near an edge which is smooth 
on the scale ro, we model the sample by a half-plane and 

set up the coordinate system so that px is the coordinate 
along the boundary and py is positive inside the half- 
plane. Then the boundary condition (54) takes the form 
z'p {px, +0) = 0. If we place the center of the saddle-point 
solution Zs shown in Fig. 4 on the edge of the sample, the 
resulting function 



Ze{p.,Py)=Zs{[{p.-p^°^f+pl]'/'^ 



(76) 



automatically satisfies not only equation (53) but also 

the boundary condition. Therefore, the expression (76) 
gives the saddle point density for the half-plane. 



One can argue that there are no other saddle point so- 
lutions for edge switching. Indeed, suppose that we have 
a solution z{px,Py) of Eq. (53) for a half-plane. Then 
we can define function z{px,py) in the entire plane, so 
that z = z{px,Py) for py > 0, and z = z{px,—py) for 
py < 0. By construction z{px,Py) satisfies Eq. (53) at 
py ^ 0. However, this procedure does not guarantee that 
the derivative z'p^ is continuous at = 0; as a result 
d'^z/dpy may have a delta- function contribution. More 
specifically, z{px,Py) satisfies the equation 



-W^z + z- 



2~z'^^{px,+mpy) 



(77) 



If in addition z{px,Py) satisfies the boundary condi- 
tion (54), i.e., Zp^{px,+0) = 0, equation (77) coincides 
with Eq. (53) everywhere in the plane. Then by construc- 
tion z{px, Py) — Zs{p), and therefore z{px,Py) is given by 
a half of the saddle-point solution Zs shown in Fig. 4 
with its center on the boundary of the half-plane. Thus, 
there are no saddle point solutions for edge switching ex- 
cept (76). 

The main exponential dependence of the mean switch- 
ing time T is given by e^^^'^. In the definition (34) of 
F{z} the integral is taken over the area of the sample. 
In the case of switching far from the boundaries it is over 
an entire plane, while for the edge switching this integral 
is over a half-plane. Therefore F is reduced by a factor of 
2 compared to the case of switching far from the bound- 
aries. Thus, instead of Eq. (56), the expression for r at 
the edge takes the form 



Te = Tg exp 



4C 



rjaiVth - V) 



1 



(78) 



The evaluation of the prefactor t* is similar to the 
one for the switching in the middle of a large sample. 
Sec. VAl. In that case we found two types of modes 
for the azimuthal part of the eigenfunctions (f)nm{p) of 
Eq. (60), namely, simrup and cos rrup. At the edge only 
the eigcnfimctions proportional to cos mip arc consistent 
with the boundary condition z'p^{px, +0) = on the di- 
mensionless density z. In the notations of Sec. VAl 
these modes correspond to m = 0, 1, 2, . . . . 

The functional F{z} is invariant with respect to the 
shifts of density z{px,Py) along the edge of the sample. 
Thus F{z} has a single zero mode .tqi; it corresponds to 
the eigenfunction with the azimuthal part cos tp. Inte- 
gration over the zero mode, in analogy with Eq. (63), is 
performed as 



(79) 



where V is the perimeter of the sample. 

To evaluate the prefactor we again use formula (36). 
Expression (41) for the numerator and the formulas for 

Dij and df/dxQo in large samples arc still applicable, 
as they were obtained in a way independent of the exact 
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form of the saddle-point density. Following the procedure 
of Sec. V A 1, we find the pref actor r* in the form 



^2C|Aoo|67i/4a5/4p 



" K 



n,m>0 



(80) 



c.f. Eq. (67). The definition of T" assumes that the 
factors corresponding to the two lowest eigenvalues, Aqo 
and Aoi, are excluded. 

In the product T" the quantum number m changes 
from to oo, while in T'" the same product is from — oo 
to oo. Thus using the fact that Xnm and Xnm are even 
functions of m, we obtain 



lnT"~ -lnT"'~ -^Inf^) 
2 IT \l J 



(81) 



see Eq. (71). 

Similarly to Eq. (72) we find the prefactor of r for the 
edge switching 



1 



C/tt 



(82) 



Note that due to the ultraviolet divergence of T" we can 
evaluate r* only up to an undetermined constant. 

Performing the same renormalization (73) of the 
threshold voltage as in Sec. VA2. one can eliminate the 
explicit dependence l^/'" of the prefactor on the mean free 
path and obtain the following expression for the mean 
switching time at the edge 



1 



exp 



7 



(83) 



The exponent in Eq. (83) is factor of 2 smaller than the 
exponent of r for the switching far from the boundaries, 
Eq. (74). Far from the threshold the exponential factor is 
dominant, and therefore edge switching is more efficient. 
To determine which switching mechanism is more effi- 
cient near the threshold, one needs to take into account 
the dependences of the prefactors in Eqs. (74) and (83) 
on the dimensions of the device. 



2. Nucleation in a comer 

In Sec. VB 1 we considered the processes of switching 
initiated near a smooth edge of the sample. In samples 
with pronounced corners, such as the devices of square or 
triangular shape, there is also a possibility of nucleation 
in a corner. As we will show, such processes may be 
more efficient than the nucleation in the interior and at 
the edges of the sample. 

We consider a corner of angle 9 < n. Similarly to the 
discussion in the beginning of Sec. VB 1, one can show 
that the saddle-point solution Zs{p) centered at the cor- 
ner both solves the equation (53) and satisfies the bound- 
ary condition (54). 



The subsequent consideration is similar to the one 
for the switching at a smooth edge of a large sample, 
Sec. VB 1. At < 77 the functional F{z} does not pos- 
sess translational symmetry with respect to the shifts of 
z{p), and therefore there are no zero modes. Due to 
the boundary condition (54) the allowed modes of the 
azimuthal part of the eigenfunction (pnmip) of Eq. (60) 
are cos{Trm(p / 9) . Then instead of Eqs. (78) and (80), we 
obtain 



27rT' 



|Aoo|^V7" 



■ exp 



4:9Cr]a 



777 




(84) 

Here A„„i are the eigenvalues of Eq. (60) with 
the boundary conditions (54), which take the form 
{d(pnm/d(p)\ip=o,e = for the corner switching. Unlike 
in Sees. V A 1 and V B 1 , here at n 1 the eigenval- 
ues Xnm are given by Eqs. (68), (F4) with m replaced by 
■Km/ 9. The product T' excludes the factor corresponding 
to the negative eigenvalue Aqo- 

Following closely the calculations of Sees. V A 1 
and VB 1, one can find the prefactor of Tc, and the ex- 
pression for the mean switching time takes the form 



1 



b^/ia \roJ 



I Y^/'^' 

— exp 



A9Qria 



777 



(85) 



One might expect that at 9 ^ it this result should 
coincide with Eqs. (78), (82) describing the edge switch- 
ing. On the other hand, the prefactors for the edge and 
corner switching are qualitatively different, since the lat- 
ter does not depend on the perimeter V. This is due to 
the fact that at < tt there is no zero mode, i.e., all 
Xnm except Aoo are positive. At ^ tt the eigenvalue 
Aoi 0, which corresponds to the appearance of a zero 
mode. In this case one needs to apply the same proce- 
dure as in Sec. VB 1, which will lead to the result (82) 
for the prefactor. 

Performing the same renormalization (73) of Vth as in 
Sees. V A 2 and V B 1, we find the expression for the mean 
switching time at a corner of angle 6 



1 



exp 



777 



(86) 



Note that ai 9 < tt the exponent of t for the corner 
switching is smaller than that for both interior and edge 
switching. This makes corner switching more efficient far 
from the threshold. 



VI. DISCUSSION 

In preceding sections we studied the mean time r of 
switching from the metastable to the stable current state 

in double-barrier resonant-tunneling structures. We cal- 
culated both the exponentials and prefactors of r for 
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switching in the small sample regime [Eq. (18)] and for 
the interior, edge, and corner switching in the large sam- 
ple regime [Eqs. (74), (83), and (86), respectively]. In this 
section we discuss the dependence of the mean switch- 
ing time on voltage for different structural parameters of 
DBRTS. 

We concentrate on the case of round samples, such as 
the ones used in the recent experiments. As we have 
shown, when the voltage V is tuned close to the thresh- 
old, the size of critical nucleus tq is large compared to 
the radius of the sample L, and the device is in the small 
sample regime. If the voltage is far from Vih, the device 
is in the large sample regime, L ^ tq. In a typical exper- 
iment T is measured in a single device for different values 
of bias. We will therefore assume that all structural pa- 
rameters and the size of the sample are fixed, and discuss 
the switching time as a function of voltage. For compar- 
ison with experiment we will not distinguish between m 
and To in this section, since the logarithmic in voltage 
corrections due to the rcnormalization of the threshold 
voltage Vth are more challenging to observe. 

Our approach is valid as long as the exponents in the 
expressions for the switching time, Eqs. (18), (56), and 
(78) are much greater than unity. To check when these 
conditions are satisfied, it is convenient to write the ex- 
ponent in Eq. (56) as 



7 



(87) 



Here we introduced a new characteristic length scale 

.3\l/4 




V 7^ 



(88) 



FIG. 5: Schematic dependence of the logarithm of the mean 
switching time r on voltage, (a) In round samples at L ^ d, 
close to the threshold there is a region of 3/2-power law de- 
pendence of Inr, followed by the region of linear dependence 
corresponding to the switching at the edge in the regime of 
large sample. In samples with pronounced corners the re- 
gion of linear dependence corresponds to the switching at the 
sharpest corner (dashed line), (b) At L ^ d, two regions of 
different linear behavior corresponding to the switching in the 
interior and at the edge of large round sample are present. In 
samples with pronounced corners these two regions are fol- 
lowed by an additional region of linear voltage dependence 
corresponding to the switching at the corner of smallest angle 
6 shown by dashed line. The slope of this linear dependence 
is smaller by a factor of n/d than that of edge switching. 



and applied the definition of given by Eq. (31b). Note 
that the length scale d depends on structural parameters 
of the device, but not on the sample size or bias. 

Similarly, the exponent of the switching time (18) in a 
small sample can be expressed in terms of d and ro as 



4ga3/2 
37V2 



6Cr2 



(89) 



where we used the fact that in round samples S = ttL^. 
This exponent is much greater than unity at ro <C 
{Ld'^ Y''^. On the other hand, the regime of small sample 
is defined by the condition tq L. Therefore, it exists 
only in sufficiently small samples, L <^ d. In this case 
close to the threshold there is a region of 3/2-power law 
behavior (18). As voltage tuned further away from Vth 
(i.e., at ro <C L), it crosses over to the region of linear 
voltage dependence of In r for the regime of large sample, 
see solid line in Fig. 5a. 

In large round samples the mean switching time t is 



given by r 



-l-Tg ^. Therefore, to find the slope of 



linear segment of the curve in Fig. 5a, one has to compare 
the rates of switching in the interior and at the edge. 



Using Eqs. (74) and (83), the ratio of the rates can be 
expressed as 



d 

-exp 



31n^ 



(90) 



A.t L ^d, this result shows that the switching always oc- 
curs at the edge rather than in the interior of the sample. 

To summarize, we found that in samples of radius L <C 
d starting at voltage difference {Vth — V) corresponding 
to ro ^ [Ld^Y/^, one first observes the region of 3/2- 
power law dependence (18) of Inr. Then, as {Vth — V) 
increases, follows the region of linear dependence (83) 
corresponding to the switching at the edge in the regime 
of large sample, see solid line in Fig. 5a. 

At 1/ ^ d the system is never in the small sample 
regime. In this case the dependence of Inr on voltage is 
linear, but it may be due to either interior or edge switch- 
ing. According to Eq. (90), at ro ^ d and very large L 
interior switching dominates. At very small ro the expo- 
nential in Eq. (90) becomes very large, and therefore the 
switching takes place at the edge. The crossover voltage 
V* between these two regions of linear dependence can 
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be determined from the condition t^^ = t- ^ applied to 
Eq. (90), 



(91) 



Thus, in these large samples the interior switching (74) 
dominates between V* and Vth, whereas at voltage below 
V* the edge switching (83) prevails. The dependence 
of Inr on voltage fov L ^ d is shown schematically in 
Fig. 5b by solid line. 

If the sample size is of order d, the dependence of In r 
on voltage can be obtained from the dependences shown 
in Fig. 5a and Fig. 5b. At L ^ d the region of 3/2-powcr 
law depedence in Fig. 5a and the interior switching region 
in Fig. 5b disappear. Thus, at L one can only observe 
the region of linear voltage dependence corresponding to 
the edge switching. 

In samples with pronounced corners the dependence 
of Inr on voltage is different due to the possibility of 
corner switching. The mean switching time r in these 
samples is given by = +t~^ +Tc^- At L ^ d an 
analysis similar to that for round samples shows that the 
region of linear dependence corresponds to the switching 
at the sharpest corner, Eq. (86). This dependence is 
illustrated by dashed line in Fig. 5a. At L ^ d, the two 
regions of interior and edge switching are followed by an 
additional region corresponding to the switching at the 
sharpest corner as {Vth — V) becomes large, see dashed 
line in Fig. 5b. 

Depending on the ratio of L and d two qualitatively 
different voltage dependences of In r are expected. To see 
whether it is possible to observe them experimentally, we 
make a crude estimate of the parameter d. Substituting 
7] = a/cb into Eq. (88) and using the estimates of 7 and 
B found in Appendix A, we get 



TrJ 
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(92) 



To obtain this expression the capacitance of the device 
per unit area was estimated as c e^n/Ep, and the en- 
ergy of the level in the well £^0 was assumed to be of 
the order of Ep. The electron density n in the well is 
typically of the order of 2 • 10^^ cm~^. The transmission 
coefficients Tj, ^ of the left and right barriers can be var- 
ied in the range from 1 to 10"**, whereas the conductivity 
measured in units of e^/fi- varies from 1 to 100. Assum- 
ing Tl ~ Tji, the low bound d ~ 20 nm is achieved at 
a ~ e^/h and Tl,r ~ 1. The upper bound d ~ 600 /zm 
is achieved by substituting the maximum value of the 
conductivity and the minimum value of the transmission 
coefficient. These estimates show that both the cases of 
L <^ d and L ^ d are experimentally achievable in mod- 
ern DBRTS, as the sample sizes range from 1 to lO^fim. 

The available experimental data^° confirm that the de- 
pendence of the mean switching time on voltage is indeed 
exponential. Based on Eq. (92) we estimate d ^ 10 /im, 
which is somewhat smaller than the radius of the sample 



L = 60 /im. Thus, one should expect the logarithm of the 
mean switching time to behave as shown in Fig. 5b. (The 
switching time r is referred to as the relocation time in 
Ref. 10.) On the other hand, it was observed in Ref. 10 
that In T bends upwards, which suggests that L <C d, see 
Fig. 5a. One of the possible explanations can be that this 
experiment was performed in superlattices, rather than 
in DBRTS studied in this paper, which makes our esti- 
mate of d unreliable. To test our theory in more detail, it 
would be interesting to carry out similar measurements 
of r in several samples of different size but with the same 
structural parameters. This will ensure that both depen- 
dences depicted schematically in Fig. 5a {L <C d) and 
Fig. 5b (L ^ d) will be observed. In addition, the ex- 
ponential dependence in Ref. 10 is not very pronounced, 
since r varies by only one order of magnitude. This sug- 
gests that T was measured rather close to the threshold, 
and therefore the data captures only the initial part of 
either linear dependence for interior switching [Fig. 5b] 
or 3/2-power law dependence. Fig. 5a. To observe the 
entire bias dependence shown in Fig. 5a or Fig. 5b, a 
measurement of r in a wider range of voltage is needed. 
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APPENDIX A: CALCULATION OF 
COEFFICIENTS A AND B IN EQ. (12) 

In this appendix we find the functions A{N) and B{N) 
in Eq. (12) in the vicinity of the threshold. We will as- 
sume that the parameter (15) is small, A <C 1. It will be 
convenient here to consider A and B as functions of the 
electron density n rather than N — nS. 

Let us write the expression for A = Jl — Jr near the 
threshold. On the upper branch of the I-V curve, the 
level in the well lies within the conduction band in the left 
lead; from Eq. (9) we obtain eV/2 + e'^N/2C + E^ > eV. 
On the lower branch, the level is below the bottom of 
the conduction band in the left lead eV/2 + Eq < eV, 
so that no current can flow through the well and N = 0. 
Therefore, in the bistable region 



eV 



e^N eV 
~2C ^ ~Y 



Eo < eV. 



(Al) 



At A < 1, it follows from Eq. (8) that e^N/2C is small 
in comparison with Ep and Eq. One can then see from 
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Eq. (Al) that ~ eV/2, and Eq. (9) results in ecj) ~ 
eV/2. Then from Eq. (10b) we find Tr = V^gnEo- The 
expression (7) for the rate Jl can also be simplified. For 
A ^ 1 and ^ gn, to first order in A the expression 
in the square brackets of Eq. (7) is {Sm/TTh'^)Ep. Using 
Eqs. (5), (7) and (10a) with all the above simplifications, 
close to the threshold A = Jl — Jr can be approximated 



as 



A{n,V) = ^1 



eV e^nA Sm 



V2gREo 



nS. 



(A2) 



The density n on the metastable and unstable branches 
of the I-V curve is found by solving the equation 
A{n, V) = 0, which reduces to the quadratic equation 

2 



'2c 



- A 



A^ 

gR 



gR 



2Eo 



El 
2Eo 



2 

en 



2c 



En 

2 ' 



0. 



(A3) 



At the threshold the two solutions for n coincide. This 
condition enables us to find the threshold voltage and 
density 



2Eo 

e 

x^/gL 
2 [gR 



gl El 



1 + — ^ 



9r El 
z El 



En 



(A4) 



(A5) 



Using Eqs. (16), (A5) and the fact that near the thresh- 
old Vr = \/2gREQ, we find the value of B{n) at the 
threshold: 



B{n,H) = ^)^\^^^El. 



(A6) 



Using Eqs. (A4), (A5) we expand A{n,V) given by 
Eq. (A2) in Taylor series near {nth,Vth) up to the first 
non- vanishing terms in n — nth and V — Vth, respectively. 
At the threshold A{n) = has only one solution, i.e., 
the first derivative with respect to n equals to zero, and 
therefore we need to expand up to the second order in n. 
The result can be presented as 



A{n) 



Bjnth) 
2 

A fan 

\9L 
A^ W 



[-a + 7(n - nthY 



El 



E'^ ' 



(A7a) 
(A7b) 



(A7c) 



Since dU/dN = —2A/B, the coefiicients a and 7 coincide 
with those used in Eq. (14). 

Assuming rectangular potential profile in the well, pa- 
rameters gL,R can be estimated in terms of the transmis- 
sion coefficients of the barriers as gL,R = Tl^r/tt. 



APPENDIX B: DERIVATION OF EQ. (21) FROM 
THE DETAILED BALANCE PRINCIPLE 



Let us consider two very close to each other points ri 
and r2 in the well. The system is in a local equilibrium, 
and the electron distributions arc given by Fermi func- 
tions. We assume electrons to be sufficiently well coupled 
to the lattice, so that the temperature T is the same ev- 
erywhere in the quantum well. Then the probabilities of 
diffusion between these two points are given by 

W^At(ri,r2;n) = ^W^,,/,(l - /,•) = ^ W^,,./,/,.e^, 

ij ij 
ij ij 

(Bl) 

Here i and j label the energy levels at positions ri and 
Y2, respectively; arc the Fermi functions, and Wij 
is the probability of transition from occupied level i to 
unoccupied level j. 

In equilibrium the transition rates satisiy the detailed 
balance condition: 



-e./T 



(B2) 



Our system is away from equilibrium, since the electro- 
chemical potential /i(r) varies with the electron density 
n(r). However, expression (B2) is still applicable for 
the relevant electron scattering processes. For example, 
in the case of elastic scattering by impurities = e^, 
and Wij = Wji due to time reversal symmetry, so that 
Eq. (B2) holds. Furthermore, one can easily check that 
for clcctron-phonon scattering expression (B2) is also 
valid, because the phonons are not sensitive to the change 
in electrochemical potential. 

Strictly speaking in the presence of electron-electron 
scattering expression (B2) is incorrect. If electron during 
the transition from state i to j scatters off an electron 
at position , the latter moves to position r'2 . Then one 
finds an additional factor of exp[{fi[ — ^2)/'?"] in the right- 
hand side of Eq. (B2). However, because the electron- 
electron interaction is screened, the distance — r'^ is 
of the order of the screening length in the well. The 
change of /i at such short distances is small compared to 
the temperature, and thus Eq. (B2) is still approximately 
correct. 

Applying expression (B2) to Eqs. (Bl) we obtain 
Eq. (21). Since during a short time interval At an elec- 
tron can only difi'use over a short distance, the above 
proof is sufficient for the purposes of Sec. HI A. 

As an additional remark, let us show that the expres- 
sion (21) also holds at larger distances. Wc consider the 
probability density Wt(rj,r/;n) of difi^usion from point 
Tj to a relatively distant point r/. Let us divide the time 
interval t into N small intervals At = t/N. Then Wt can 
be represented in terms of VFa* in the following way: 

N 

Wt(rj,r/;n) = ]J WAt{r-k,rk+i;n) dvk+i, 
k=i •' 
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where ri = and rpf+i = Vf. The distances between the 
points Tfe and r^^i are small, so that the expression (21) 
is applicable. Since at small distances A/i <C T, we can 
expand Eq. (21) up to the linear terms in A^/T. Using 
this expansion we can rewrite each integrand in the above 
expression in terms of Wai (ffe+i , r/c ; n) . Then evaluating 
the product over k we obtain Eq. (21). This completes 
the prove. 



Changing the variables to y = yi — y2 and f = (j/i + 
2/2)72, and using the fact that jdv = Ly, we obtain 



4TAt 



Jdvr^WAt{\r\;n). (C4) 



Finally, comparing Eqs. (24) and (C4) we find G = 



APPENDIX C: CALCULATION OF CONSTANT 
G IN EQ. (24) 

In this appendix we find the constant G in Eq. (24) for 
an arbitrary scattering mechanism. This is accomplished 
by expressing G in terms of conductivity a. 

If a small electrochemical potential gradient is applied 
in the x-direction, it gives rise to an electric current 



(CI) 



where Ly is the width of the sample. 

Let us find the expression for the current along x-axis 
at j; = in terms of the transition probability density 
Wa*- It is given by the difference in the number of elec- 
trons crossing the line .x = from left to right and in the 
opposite direction in unit time, namely, 

dy2 I dxi dx2 
J-00 Jo 

X [WAt (ri , r2 ; n) - WAt (r2 , r 1 ; n)] . (C2) 



In equilibrium, i.e. at dji/dx = 0, the expression for the 

difference of probability densities in the second line of 
Eq. (C2) vanishes. Away from equilibrium it can found 
by using the "detailed balance" expression (21): 

WAt{ri,r2;n) - VFAt(r2, ri; n) ~ ^ WAt(ri, rz; n). 

Expanding fii — fi2 ~ {xi — X2)d[i/dx, one can see that 
the linearized form of Eq. (C2) reproduces Eq. (CI) with 
the conductivity given by 



a = 



I dyi / dy2 / dxi / dx2 
Jo Jo J-00 Jo 

x{x2-xi)WAtiri,r2;n). (C3) 



LyTAt 



It is important to note that this expression is taken 
in the hmit dfi/dx = 0, so that WAt in Eq. (C3) is an 
equilibrium quantity. Therefore Wai depends only on 
the distance between ri and r2, i.e., T4^At (ri , r2 ; ?^) = 
W^At(|i"i — r2|;n). Then substituting new variables x = 
X2—X1 and u = {xi +X2)/'^ into the integral in Eq. (C3), 
and integrating over u, we find 

Ly I'Ly I'OO 

dyi I dy2 I dx 



2LyTAt J a 







xx^WAtiVx^ + (2/1 -2/2)2; 



APPENDIX D: STATIONARY SOLUTION OF 
EQ. (28) NEAR THE THRESHOLD 

In this appendix we discuss the stationary solution 
Ps{n} of the Fokker-Planck equation (28) near the 

threshold. At bias near Vui function b{n) can be ap- 
proximated by a constant b = b{nth)- Then the equation 
for Ps{n} takes the form, 



u'{n) 



5 o 2(7 o 6 

on e b on 



P,{n} = 0. (Dl) 



Here 77 = a job and 2o/6 = lAjB = -u'{n), c.f. Eq. (13). 

It is convenient to present Ps{n} in terms of a func- 
tional Fs{n}, such that Ps{n} = exp(— i^gjn}). Then 
Eq. (Dl) takes the following simple form 

2^V2/(r)+j/(r) = /(r), 



where we introduced y{r) = SFg/Sn and /(r) 
2riV'^n. 

Solution of this equation is given by 



y{r) = J dr' f{r')g{r-r'), 



(D2) 
u'{n) - 

(D3) 



where Q is presented in terms of the modified Bessel func- 
tion Kq as 



2aT 



(D4) 



At low T the characteristic size of the Green's func- 
tion Q is very small, so that Q can be approximated by a 
(5-function. Then Eq. (D3) greatly simplifies, 

— F,{n} = u'(n(r)) - 27)V^n{v). (D5) 
on 

The solution of this equation reproduces Eq. (30). 

The exact criterion for validity of Eq. (D5) is given 
by the condition r^, where tq is the characteristic 

size of the function —u'{n) + 2rf\/'^n, see Eq. (31b). Af- 
ter substitution of the parameters of the problem from 
Eqs. (A7) and (31b) this criterion takes the form: 

r< a2(i + a) -^^ 



\9rJ El'^^e{Vth-V) 



(D6) 



To estimate the right-hand side of (D6) we take the pa- 
rameters A ~ 1, E'f ~ Eq, e{Vth-V) < Ef and gL ~ Qr- 
Then the criterion (D6) reduces to T <i; Ep- Therefore, 
one can neglect the temperature term in Eq. (28) unless 
the structure is strongly asymmetrical, so that gL ^ gn- 
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APPENDIX E: PROPERTIES OF THE SADDLE 
POINT SOLUTION Zs{p) 

In this appendix wc derive several relations between 
integrals involving Zs{p). Our goal is to express the in- 
tegrals / (Vzs)'^ dp and / Zg dp in terms of ( defined by 
Eq. (57). 

Integrating Eq. (53) over the infinite plane and using 
the fact that Zs{p) decays rapidly at large p, we find 



:dp = / Zg dp. 



(El) 



To express the integral in Eq. (64) in terms of the con- 
stant C, we transform it as 



{Vzsfdp = - J ZsV^Zgdp 

= -JzUp + jzldp, (E2) 

where we used Eq. (53) to obtain the second line of (E2). 

To express one of the integrals in the second line of 
Eq. (E2) in terms of the other, we take advantage of the 
azimuthal symmetry of the saddle point solution. Multi- 
plying Eq. (55) by p^z'^ and integrating over p, we find 



r 2 d (zl zl\ , „ 



The first term in this equation vanishes, since z^ oc at 
p — !■ cx) . The second term can be simplified by integration 
by parts, resulting in 



Za dp. 



(E3) 



Using Eqs. (El) (E3) and (57) we find the following 
expressions for the integrals in Eqs. (64) and (F6), 



{Wzsfdp = 2C, 



Zs dp = 4C. 



(E4) 



APPENDIX F: SOLUTIONS OF EQ. (60) 

In this appendix we find the eigenvalues of continuous 
spectrum of the Schrodinger equation (60). We consider 
a round sample of dimensionless radius R = ^ Sj nr^ 
with the critical fluctuation situated in the center. Note 
that since we are interested in the case of large samples, 
the size of the critical fluctuation is small compared to 
the sample size, i.e., i? 1. 

The potential — 2zs(/3) is azimuthally symmetric, so it 
is convenient to solve equation (60) in polar coordinates. 
Separating the variables in as (j)nm{p) = Qnm{p)'^m{^)i 
we can write the equation for the radial part as follows 



d^ I d 

dp-' pdp p^ 



Qnm{p) — QnmQri 



(Fl) 



where q. 



2 

nm 



1. This equation is subject to two 



boundary conditions: Qnm{p) is finite at the origin and 
%™(i?) = 0. 

Let us flrst consider an infinite sample. In the absence 
of the attractive potential — 2zs, the finite at the origin 
solutions to Eq. (Fl) are the Bessel functions of the first 
kind Jm{QP)- Their asymptotic behavior at p ^ 00 is 



Jmiqp) - \ COS 

nqp 



(F2) 



In the presence of the attractive potential the asymptotic 
form of the radial part of the eigenfunction modifies as 
follows 



Qm{p) - -1/ COS 

" Tvqp 



■ (F3) 



Here 5m{q) is the scattering phase shift due to the at- 
tractive potential. 

For our purposes we only need the expression for Sm at 
large wavevectors q. At q^ 1 the phase shifts (5m <C 1, 
and can thus be found in Born approximation. 



Sm{q) = TT / Zs{p)Jl^{qp)pdp, 
Jo 



(F4) 



see also Eq. (14) in Ref. 36. Note that 6m is indeed small 
at g » 1, because oc 1/q. 

In a finite sample the wavevectors are quantized. 
Using the asymptotic form (F3) and the boundary con- 
dition Q'nmi^) = 0' fi^d 



qnm — 



6m iflnra) 



(F5) 



where fi is given by n -|- 1/4 if m is even, and by n + 3/4 
if m is odd, with n being a nonnegative integer. Then 
the eigenvalues A„„i are given by 1 + (j^„„. We use this 
result in Sec. V A 1 to calculate ■ This product 

is dominated by the factors with large qnm- Therefore in 
Eq. (68) we approximate fi by the radial quantum number 
n and the argument of 5m{qnm) by T^n/R. 

In addition, in Sec. V A 1 we need an expression for the 
sum of the phase shifts (F4) over the azimuthal quantum 
numbers m. In the right-hand side of Eq. (F4) only the 
Bessel functions Jm{qp) depend on m. Since the sum of 
J^iqp) over m equals unity,^'' we find 



,{p)pdp = 2C. 



(F6) 



We used Eq. (E4) to express the above integral in terms 
of the constant (. 

This result can also be derived by means of the Friedel 
sum rule which states that the sum of the phase shifts in 

the left-hand side of Eq. (F6) is given by ttN, where N is 
the average number of levels in the attractive potential 



21 



U = —2zs{p). Since the two-dimensional density of states Combining this expression with the Priedel sum rule we 
1/2 = 1/477, we find reproduce the result (F6). 

N = ^ (2zs{p)dp. (F7) 
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